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ABSTRACT 


This  report  describes  an  investigation  of  the  ability  to  predict  long-term 
orbital  elements  for  medium-height  (satellite  orbiting  the  earth  at  altitudes 
above  the  regions  of  appreciable  drag),  earth-orbiting  satellites.  The  partic¬ 
ular  types  of  satellites  studied  are  those  having  nearly  circular,  polar 
orbits.  The  primary  element  of  interest  is  orbital  inclination. 

A  theoretical  equation  for  predicting  the  change  in  inclination  caused  by 
the  gravitational  attraction  of  the  sun  and  moon  is  given.  Expressions  for 
the  orbital  e le ments  of  the  sun  and  moon  required  in  this  inclination- 
predicting  equation  are  developed.  A  FORTRAN  program  listing  and  compu¬ 
tation  procedure  for  determining  the  inclination  of  a  satellite  at  any  time  after 
a  chosen  epoch  is  in  the  appendix. 

Nine  satellites  were  selected  for  study.  For  these  satellites,  compari¬ 
son  plots  of  observed  and  predicted  inclinations  over  a  3-year  time  interval 
are  shown.  The  agreement  between-  predicted  and  observed  inclinations  is 
very  good  for  most  of  the  satellites.  Theoretical  plots  of  the  inclination  of 
Satellite  902  over  time  intervals  of  up  to  300  years  are  included. 


PROBLEM  STATUS 

This  is  an  interim  report;  work  is  continuing  on  the  problem. 
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INCLINATION  PERTURBATIONS  OF  POLAR-ORBITING  SATELLITES 


INTRODUCTION 

With  the  advent  of  modern  satellite  navigation  systems,  the  ability  to  generate  a 
satellite  ephemeris  which  is  valid  for  a  long  period  of  time  has  become  important.  Much 
effort  is  now  being  directed  toward  the  development  of  techniques  for  accurately  pre¬ 
dicting  the  position  of  a  satellite  at  any  instant  of  time.  Sophisticated  models  of  the  earth’s 
gravitational  field  and  of  the  upper  atmosphere  have  been  developed  to  explain  the  observed 
complex  motion  of  an  orbiting  satellite.  Much  work  is  left  to  be  done  in  this  area  of 
predicting  satellite  orbits. 

In  many  '-ases,  knowledge  of  the  short-term  satellite  position  is  not  the  only  important 
factor.  Long-term  orbit  predictions  can  be  as  important  in  applications  wherein  station 
keeping  is  a  necessity  or  wherein  knowledge  of  the  maximum  deviations  of  orbital  ele¬ 
ments  from  their  initial  values  is  required.  Examples  of  such  applications  are  satellites 
in  synchronous  or  subsynchronous  orbit  or  satellites  required  to  traverse  a  specified 
constant  or  nearly  constant  ground  track  over  a  long  period  of  time. 

In  line  with  the  above  purposes,  it  is  desirable  to  check  the  validity  of  existing  theory 
on  the  long-term  prediction  of  orbits  and  to  ascertain  the  important  factors  which  must 
be  controlled  or  compensated  for  in  order  to  maintain  any  desired  satellite  ephemeris. 


COORDINATE  SYSTEM  AND  TERMINOLOGY 

The  position  of  a  satellite  is  described  by  the  six  conventional,  elliptical,  orbital 
elements.  The  position  of  the  orbital  plane  along  the  equator  is  defined  by  the  right 
ascension  of  the  ascending  node  fl  and  is  measured  eastward  from  Aries  {the  vernal 
equinox).  The  angle  between  the  equatorial  and  orbital  planes  is  defined  as  the  orbital 
inclination  i ,  The  position  of  the  ellipse  in  its  own  plane  is  defined  by  the  argument  of 
perigee  *>.  The  size  and  obliquity  of  the  orbit  are  defined  by  the  semimajor  axis  a  and 
the  eccentricity  <•,  respectively.  The  instantaneous  position  of  the  satellite  in  the  orbital 
plane  is  defined  by  the  mean  anamoly  tf  measured  from  the  argument  of  perigee.  Figure  1 
illustrates  the  above  definitions. 


THEORY 

Of  the  six  orbital  elements  defined  in  the  last  section,  orbital  inclination  is  one  of 
the  easiest  to  work  with  from  the  standpoint  of  both  ease  of  prediction  and  physical  inter¬ 
pretation.  Therefore,  in  the  following  work,  inclination  will  be  the  element  of  interest. 

For  high,  circular-orbiting  satellites,  the  primary  forces  affecting  satellite  motion 
are  variations  in  the  gravitational  field  of  the  earth  and  the  gravitational  attraction  of  the 
sun  and  moon.  For  satellites  with  a  large  surface-area-to-mass  ratio,  solar  radiation 
pressure  may  also  be  important.  The  relative  importance  of  the  above  influences  vary 
according  to  the  orbital  element  considered.  For  example,  the  most  prominent  influence 
on  right  ascension  and  argument  of  perigee  is  related  to  irregularities  in  the  earth’s 
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Fig.  1  -  Elliptical  coordinate  system 


gravitational  Held,  while  the  major  effect  on  orbital  inclination  is  caused  by  the  lunisolar 
forces.  Perturbations  in  inclination  are  only  indirectly  caused  by  the  earth  through 
changes  in  the  other  satellite  elements  (e.g.,  right  ascension,  argument  of  perigee,  etc.). 

Theoretical  expressions  for  the  rates  of  change  of  orbital  elements  for  satellites 
under  the  influence  of  lunisolar  gravity  forces  have  been  developed  by  several  authors 
(1,2,3).  A  report  by  G.E.  Cook  (1),  for  example,  lists  expressions  for  computing  first- 
order  rates  of  change  in  orbit  parameters.  An  extension  of  Cook’s  work  in  a  paper  by 
Myrna  Lewis  (2)  gives  the  pertinent  expressions  in  a  form  more  suitable  for  computation, 
the  equation  listed  for  inclination  is  the  one  used  to  obtain  predictions  in  this  report.  One 
of  the  most  comprehensive  papers  to  date  on  perturbations  of  satellite  orbits  is  by  R.H. 
Gooding  (3).  This  report  contains  higher-order  terms,  as  well  as  the  primary  ones 
appearing  in  the  above  reports,  for  computing  the  rates  of  change  of  satellite  orbital 
elements.  Several  other  possible  perturbations  are  also  discussed  in  his  report.  Since 
the  purpose  of  this  report  is  not  to  study  the  theory  of  satellite  perturbations,  but  to 
determine  its  validity  and  usefulness,  only  the  appropriate  results  of  these  theoretical 
developments  will  be  given. 

The  rate  of  change  of  orbital  inclination  may  be  predicted,  fairly  accurately,  by 
Eq.  118  in  Ref.  2.  If  the  satellite  orbit  of  interest  may  be  regarded  as  circular,  those 
terms  in  the  equation  multiplied  by  e2  may  be  eliminated.  Also,  for  polar-orbit  satellites, 
terms  involving  cos  i  may  also  be  omitted.  The  resulting,  simplified  version  of  this 
equation  is  given  below: 

ft  -  — '^!)(sin  i){sin4  T[sin  2  (*  *  --  *  M 

-  —sm  2  i^p  t  (O',  ♦  -yj  -^L  sin  ^2 p  f  2 u(,  t 


(1) 
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(Equation  (1)  Continued) 

♦  cos4  [sin  2  (/?  - 


-d  sin  2 
2 


fa  ~  ud  ~ 

sin  ^2ft  -  2 to',  -  3VdjJ  *1  2sin2  *"rf  sin2  2ft  [l  .  3  e (  cos  V dj 


(1) 


In  the  above  equation,  mean  orbital  elements  are  to  be  used  (short-period  terms, 
i.e.,  those  related  to  time  intervals  of  less  than  one  satellite  revolution  are  integrated 
out).  Assumptions  used  in  the  development  of  this  equation  are  that  the  perturbing  body 
does  not  move  during  one  revolution  of  the  satellite,  and  that  the  semimajor  axis  of  the 
satellite  orbit  is  much  smaller  than  the  distance  from  the  earth  to  the  moon-at  least 
smaller  than  a  tenth  of  the  ear th -moon  distance.  In  this  expression,  variables  having  a 
subscript «/  refer  to  the  orbital  elements  of  the  disturbing  body  (sun  or  moon).  Two 
quantities  appearing  in  the  equation  need  some  clarification.  The  variable  n  is  the  right 
ascension  (n)  of  the  satellite  if  (h,dt  due  to  the  sun  is  being  calculated,  and  ft  n  -  0(/  if 
ill  ,lt  due  to  the  moon  is  being  computed  (ftd  for  the  sun  is  by  definition,  zero).  The  con¬ 
stant  K  is  defined  as 


K  ;  (2) 

where  <*d  is  the  product  of  the  universal  gravitation  constant  G  and  the  mass  of  the  dis¬ 
turbing  body.  The  remaining  variables  in  the  equation  require  no  special  mention  but  are 
defined  in  the  glossary  at  the  end  of  the  report. 


APPLICATION 

In  order  to  use  Eq.  (1)  for  inclination  prediction,  the  orbital  elements  of  the  sun  and 
moon  and  certain  satellite  parameters  must  be  determined.  The  elements  of  the  sun 
present  little  problem.  They  are  directly  available  from  “The  American  Ephemeris  and 
Nautical  Almanac”  (4)  for  any  particular  year  of  interest.  Approximations  of  the  equa¬ 
tions  listed  in  Ephemeris  for  predicting  orbital  elements  of  the  sun  are  given  below. 


These  equations  yield  valid  predictions  over  several  years  (from  1966). 

L  (mean  longitude)  =  279°6967  +  0?98564734  d  (3) 

r(mean  longitude  of  perigee)  =  281^22083  +  0!00004707  d  (4) 

L  -  r  (mean  anomaly)  =  358°47583  +  (L9856002  d  (5) 

es  (eccentricity)  =  0.01673  (  6) 

<  (mean  obliquity  or  inclination)  =  23^4432  (7) 


In  the  above  equations,  j  is  in  units  of  ephemeris  days  measured  from  the  epoch  1900 
January  0.5  E.T. 

Similar  expressions  exist  for  the  orbital  elements  of  the  moon  and  are  listed  here. 

C  (mean  longitude)  =  270?434  +  13tl7639653  d  (8) 

r'  (mean  longitude  of  perigee)  =  334°3294  +  0°,11140408  d  (9) 
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Q  (right  ascension)  =  259?1833  -  0?05295392  d  (10) 

cy  (eccentricity)  =  0.054900489  (11) 

Unfortunately,  these  elements  are  measured  from  the  mean  equinox  of  date  along  the 
ecliptic  and  then  along  the  lunar  orbit.  Therefore,  they  must  be  transformed  to  the  coor¬ 
dinate  system  defined  at  the  beginning  of  this  paper  (see  Fig.  1).  This  is  easily  accomp¬ 
lished  by  the  use  of  spherical  trigonometry.  Figure  2  (a  modified  form  of  a  diagram 
appearing  in  Cook’s  paper  (1))  illustrates  the  definitions  of  the  moon’s  elements  listed 
above  and  the  required  transformed  lements.  Expressions  for  obtaining  the  desired 
elements  derived  with  the  aid  of  Fig.  2  are 


i  (inclination)  =  cos'^oosr  cos  a  -  sin  <  sin  a  cos  Q  ) 

and 

» „  (argument  of  perigee)  =  o  .  r  -  0 , 

where o  = 

sin- i/sin  <  s»n  Q  \  ,  s  obliquity  of  the  ecliptic 

\  sin  i  f 

(12) 

n  (righc  ascension)  =  sin- 

■  i  (  ?.a  )  where  a  is  the  inclination  of  the  lunar 

\  sin  1  j 

(13) 

orbit  to  the  ecliptic 

L  (mean  longitude)  =  <f  -  Q  ♦  o 

(14) 

u  (mean  anomaly)  -  C  -  r 

(15) 

efc.  (eccentricity)  =  same  as  Eq.  (11). 

With  the  aid  of  all  of  the  above  expressions,  we  can  now  predict  the  orbital  elements 
of  the  sun  and  moon  for  any  desired  time.  The  only  remaining  information  that  is  needed 
is  the  orbital  elements  (other  than  inclination,  which  is  required  only  at  the  chosen  epoch) 
of  the  satellite  as  a  function  of  time.  These  elements  need  not  be  known  to  great  accuracy 
because,  with  the  exception  of  mean  motion  and  eccentricity,  each  element  appears  only 
as  a  sin-cos  argument  in  the  prediction  equation  for  inclination.  Two  of  the  satellite 
elements,  right  ascension  and  argument  of  perigee,  can  be  predicted  with  sufficient 
accuracy  (i.e.,  sufficient  for  use  in  predicting  inclination)  by  using  the  following  expres¬ 
sions  (modified  forms  of  equations  appearing  in  the  Space  Navigation  Handbook  (5)  for 
their  rates  of  change: 


n  =  2"nJn cos  ‘ 

a2(  1-e2)  2 

and 

<0 .  2”<(2-|sin2  ') 
o2(  1-e2)  2 

where  jn  is  the  nth  zonal  harmonic  of  the  earth.  For  the  accuracy  required  of  the  pre¬ 
dictions,  only  the  major  harmonic  J2  need  be  used.  If  n  is  in  rev/day  and  a  in  units  of 
earth  radii,  then  the  units  of  the  results  are  rad/day.  Given  the  right  ascension  and 
argument  of  perigee  at  any  epoch,  their  values  on  succeeding  days  can  be  determined  by 
simply  adding  the  cumulative  rates  of  change  to  the  initial  value.  The  eccentricity  and 
the  semimajor  axis  values  (also  mean  motion,  which  may  be  derived  from  the  semimajor 
axis),  if  changing  appreciably  during  the  required  prediction  time  interval,  must  also  be 
updated.  It  should  be  noted  that  the  mean  anomaly  of  the  satellite  is  not  required,  as 


(16) 

(17) 
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Fig.  2  -  The  orbital  elements  of  the  moon 
(modified  from  Ref.  1) 


it  has  been  eliminated  through  integration  in  the  development  of  Eq.  (1).  Equation  (1)  gives 
the  average  rate  of  change  of  inclination  over  one  complete  orbit. 


PREDICTION  PROGRAM 

A  FORTRAN  program  for  computing  the  orbital  inclination  of  a  satellite  v/as  written 
incorporating  the  previous  equations.  This  program  was  developed  with  the  intent  of 
predicting  inclination  only  for  those  satellites  in  high,  circular,  polar  orbits;  thus,  all  of 
the  terms  in  some  of  the  previous  equations,  particularly  Eq.  (1),  are  not  necessary. 

The  general  procedure  in  this  program  is  first  to  undate  the  elements  of  the  sun,  and 
then  to  calculate  the  rate  of  change  of  inclination  due  to  the  sun.  Next,  the  moon’s  ele¬ 
ments  are  updated,  and  the  rate  of  change  of  inclination  due  to  the  moon  calculated.  The 
net  rate  of  change  is  then  the  linear  combination  of  the  two  calculated  rates  of  change. 
Finally,  the  satellite  elements  are  updated  in  preparation  for  the  computation  of  inclina¬ 
tion  for  the  next  unit  of  time.  This  process  is  then  repeated  until  the  desired  time  span 
of  predictions  is  reached. 

The  input  data  required  by  the  program  includes  the  initial  satellite  orbital  elements 
with  all  angles  in  radians,  semimajor  axis  in  earth  radii,  and  mean  motion  in  deg/day 
(as  noted  previously,  the  mean  anomaly  of  the  satellite  is  not  required).  The  number  of 
ephemeris  days  from  1900  Jan  0.5  E.T.  to  the  desired  starting  point  in  time  is  also 
required.  Another  time  factor,  which  has  not  been  previously  discussed,  appears  in  the 
program.  The  variable  T,  which  is  d/36525,  is  the  time  from  1900  Jan  0.5  E.T.  to  the 
prediction  epoch  in  Julian  centuries.  It  is  only  required  for  prediction  over  a  time  interval 
of  many  years.  It  appears  in  the  more  complete  expressions  for  the  orbital  elements 
of  the  sun  and  moon  in  Ref.  4. 

The  time  interval  between  successive  predictions  and  corrections  to  all  orbital 
elements  is  set  for  1  day  in  the  program.  This  can  easily  be  changed,  however,  bv  simply 
modifying  the  factors  added  to  </  and  T  in  the  prediction  program.  The  second  index  of 
the  DO  loop  may  also  be  altered  to  allow  coverage  of  any  desired  prediction  time  span. 

The  appendix  contains  the  FORTRAN  listing  of  this  prediction  program. 
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PREDICTIONS 

It  was  decided  that  the  best  approach  to  evaluating  the  prediction  capabilities  of  the 
program  was  to  compare  computed  inclinations  on  certain  existing  satellites  with  avail¬ 
able  data  obtained  from  the  U.S.  Naval  Space  Surveillance  System.  To  examine  this 
capability,  nine  polar-orDiting  satellites  were  chosen.  The  epoch  for  predictions  and 
data  was  chosen  to  be  January  1,  1966  E.T.  This  epoch  was  chosen  because  the  data 
prior  to  this  time  was  not  of  sufficient  accuracy  to  be  useful  for  accomplishing  the  stated 
task.  Improvements  in  the  WS-434,  the  U.S.  Naval  Space  Surveillance  System  in  late  1965, 
vastly  improved  the  data.  The  mean  orbital  elements  of  the  nine  selected  satellites  on 
January  1,  1986  appear  in  Table  1. 


Table  1 

Orbital  Elements  of  Nine  Satellites  on  January  1,  1966 


Satellite 

Period 

(min) 

Period 

Rate 

(min/day) 

Inclination 

(deg) 

Eccentricity 

Semimajor 

Axis 

(earth  radii) 

Right 

Ascension 

(deg) 

Right 

Ascension 

Rate 

(deg/day) 

Argument 
of  Perigee 
(deg) 

Argument  of 
Perigee  Rate 
(deg 'day) 

509 

99.0553 

-0.00005 

90.658 

0.00279 

. 

1.1114 

f" . . 

253.602 

C.079 

130.206 

-3.427 

671 

107.3350 

-0.00000 

89.900 

0.00441 

1.1725 

4.956 

-0.010 

85.880 

-2.844 

704 

107.1210 

-0.00000 

89.947 

0.00331 

1.1710 

98.706 

-0.005 

290.443 

-2.857 

801 

103.0663 

-0.00001 

90.505 

0.00651 

1.1412 

42.382 

0.055 

102.1S2 

-3.125 

902 

106.6006 

-0.00000 

89.900 

0.00168 

1.1672 

327.698 

-0.010 

151.513 

-2.890 

959 

108.3125 

-0.00000 

89.969 

0.00396 

1.1651 

143.377 

-0.003 

267.878 

-2.908 

1314 

111.5066 

-0.00000 

90.205 

0.00279 

1.2027 

218.570 

0.019 

242.541 

-2.602 

1420 

106.9138 

-0.00091 

89.976 

0.00788 

1.1695 

311.469 

-0.002 

114.996 

-2.871 

1514 

108.1308 

-0.00000 

90.012 

0.00694 

1.1783 

355.603 

0.001 

267.997 

-2.796 

Several  prediction  runs  were  made  for  the  nine  satellites  when  the  program  debugging 
was  completed.  One  of  the  first  plots  made  of  theoretical  and  observed  data  for  Satellite 
902  (the  10-kg,  14 -in. -diameter  NRL  calsphere  (6)  is  shown  in  Fig.  3.  The  observed 
and  predicted  data  appear  to  agree  well  during  most  of  1966  and  early  1967,  but  the  two 
curves  tend  to  drift  apart  in  late  1967  and  continue  to  separate  during  1968.  The  initial 
conclusions  made  on  the  basis  of  this  plot  were  that  either  an  error  existed  in  the  com¬ 
puter  program  or  some  factor  affecting  inclination  either  had  not  been  incorporated  into 
the  prediction  theory,  or  had  been  incorrectly  included. 

After  much  labor  in  checking  the  computer  program  for  errors,  and  checking  the 
derivations  of  both  the  theoretical  equations  developed  by  Myrna  Lewis  (2)  and  the  expres¬ 
sions  for  predicting  the  sun,  moon,  and  satellite  elements,  it  was  concluded  that  no  errors 
existed  in  the  theory  being  used  or  in  its  application.  Other  factors  not  included  in  the 
theory  which  may  affect  inclination  were  also  considered  and  found  to  be  insignificant 
for  the  type  of  errors  found  (further  mention  of  these  neglected  factors  will  occur  later). 
The  only  other  item  left  to  check  was  the  WS-434  data. 

After  studying  available  literature  on  WS-434  computation  methods  and  through 
communication  with  presonnel  at  Space  Surveillance  System  headquarters  at  Dahlgren, 
Virginia,  it  was  finally  determined  that  the  problem  was  with  the  coordinate  system  ref¬ 
erences  used.  The  theoretical  data  were  referenced  to  the  true  equinox  of  dale,  while  the 
WS-434  data  were  referenced  to  two  different  mean  equinoxes  during  the  3-year  time 
interval  considered.  The  WS-434  data  collected  during  1966  were  referenced  to  the  mean 
equinox  of  January  1,  1966.  All  data  collected  since  the  end  of  1966  have  been  referenced 
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Fig.  3  -  Theoretical  and  observed  inclination  for  Satellite  902  before 
correction  for  precession  and  nutation 


to  the  mean  equinox  of  January  1,  1967.  Since  a  true  equinox  reference  is  required,  it 
was  necessary  to  correct  the  data  by  adding  the  nutaticnal  effects  and  the  cumulative 
processional  effects  occurring  after  each  reference  epoch. 


ADJUSTING  DATA  COORDINATE  SYSTEM 

General  precession  consists  of  motion  of  the  vernal  equinox  backwards  along  the 
ecliptic  at  the  approximate  rate  of  507  3  per  year  (the  full  period  is  about  26,000  years). 
This  is  the  only  type  of  precession  of  any  significance  for  the  purposes  of  this  paper 
(there  exist  other  precessional  terms  of  much  smaller  magnitude).  There  are  also  two 
nutational  terms  affecting  the  position  of  the  reference  or  the  orbital  elements  of  the 
satellite,  sun,  and  moon.  Nutation  in  longitude  involves  motion  of  the  equinox  along  the 
ecliptic,  as  does  precession.  Nutation  in  obliquity  causes  a  change  in  the  obliquity  of 
the  ecliptic  (inclination  of  the  sun’s  apparent  orbit).  Precessional  and  nutational  effects 
are  caused  by  the  gravitational  attraction  of  the  sun  and  moon  on  the  equatorial  bulge  of 
the  earth. 
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The  motion  of  the  coordinate  system  reference  and  its  effect  on  satellite  inclination 
is  illustrated  by  Fig.  4  for  the  case  of  a  satellite  having  a  right  ascension  in  the  first 
quadrant.  With  the  aid  of  this  figure,  an  expression  may  be  developed  for  calculating  the 
change  in  inclination  M,  for  any  particular  precession  or  nutation,  which  must  be  added 
to  the  observed  inclination  data.  From  a  computational  point  of  view,  it  is  not  necessary 
to  have  only  a  single  equation  for  calculating  Af.  A  set  of  short,  ordered  equations  is  of 
equal  value,  and  is  much  easier  to  develop.  Such  a  series  of  equations  for  the  case 
illustrated  in  Fig.  4  follows: 

cos  0  ■=  cos  (  cos  i  f  sio  t  sin  i  cos  fl,  (18) 

sin  0  =  sin  (arccos  (cos  6)),  (19) 


and 


sin  R  =  1  si.'- 1*.  (20) 

sin  0 

Re  arcsin  (sin  R),  (21) 

cos  (i  <■  A i2)  =  cos  ({  *  AO  cos  0  -  sin  ( i  +  Af )  sin  0  cos  (R  ♦  P),  (22) 

i  +  Aij  =  arccos  cos  (i  +  A / 2) .  (23) 

Ai2=  Ai  =  i.  Aij  -  i.  (24) 


Care  must  be  taken  when  using  this  method  of  analysis  in  connection  with  spherical  trig¬ 
onometry,  as  the  same  equations  may  not  be  valid  for  angles  in  quadrants  other  than  the 
one  being  considered. 


Fig.  4  -  Precession  and  nutation  effects 
on  inclination 
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A  simple  equation  is  given  in  a  paper  by  Newton  (7)  for  determining  the  apparent 
change  in  inclination  due  to  precession  and  nutation.  It  is  listed  here  for  reference: 


\/  -  -  (P  •  S'l)  sin  /,  5in  0S  ♦  (\/s>  cos  (ljt 


(25) 


where 


P  =  total  precession  since  an  epoch, 
NL  =  nutation  in  longitude, 
ls  =  obliquity  of  the  ecliptic, 
ns  =  right  ascension  of  the  satellite, 


and 


Ms  =  nutation  in  obliquity. 

This  equation  has  been  evaluated  experimentally  by  comparison  with  inclination  correc¬ 
tions  obtained  using  the  set  of  equations  given  above,  and  it  appears  to  be  a  very  good 
approximation. 

The  values  of  the  nutational  terms  may  be  found  for  any  given  day  in  “The  American 
Ephemeris  and  Nautical  Almanac”  (4).  It  is  much  simpler  to  calculate  them,  however, 
when  their  values  are  required  in  a  computer  program.  Equations  for  calculating  the 
approximate  values  of  nutation  in  longitude  and  obliquity  are  given  by  Smart  (Ref.  8,  p.  247). 
The  major  terms  of  these  expressions  are 


It  H  »»  ti 


NL  =  -17.234  sin  Q  +  0.209  sin  (2Q)  -  1.272  sin  (2b)  -  0.204  sin  (2C) 

(26) 

and 

V  =  9.210  cos  Q  +  0.551  cos  (2b)  +  0.088  cos  2€, 

(27) 

where 

NL 

is  nutation  in  longitude 

\< 

is  nutation  in  obliquity 

L 

is  the  sun’s  mean  longitude. 

i 

is  the  moon’s  longitude. 

and 

0 

is  the  moon’s  right  ascension. 

Figure  5  shows  a  plot  of  the  theoretical  inclination  and  the  observed  inclination 
corrected  using  the  scheme  developed  above  for  Satellite  902.  The  agreement  between 
two  curves  is  now  excellent  over  the  3-year  time  interval  considered.  All  significant 
perturbing  influences  seem  to  have  been  incorporated  into  the  prediction  theory.  Figures  6 
through  13  show  similar  plots  for  the  eight  other  satellites  studied.  Most  of  the  satellite 
predictions  hold  as  well  as  they  do  for  Satellite  902.  Three  of  the  satellites  (numbers  704. 
959,  and  1314)  exhibited  a  peculiar  behavior  near  the  end  of  1966  and  in  early  1967.  During 
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this  period,  the  predicted  and  observed  data  (corrected  for  precession  and  nutation)  plots 
tend  to  separate  and  then  remain  nearly  parallel  during  the  remainder  of  the  3 -year  period. 
The  exact  cause  of  this  effect  is  not  clear;  however,  the  peculiar  behavior  seems  to  be 
associated  with  the  WS-434  orbital  elements  for  these  satellites  rather  than  with  the 
theory.  There  are  several  reasons  for  this  conclusion.  First,  since  the  two  curves 
(theoretical  and  observed)  remain  nearly  parallel  after  the  beginning  of  1B67,  simply  chang¬ 
ing  the  epoch  to  this  data  would  allow  nearly  perfect  prediction  for  the  remaining  two 
years.  Also,  the  three  satellites  exhibiting  this  effect  should  be  well  behaved  In  their 
respective  orbits.  No  known  active  propulsion  sources  exist  for  these  satellites,  and 
atmospheric  drag,  solar  radiation,  pressure,  and  resonance  effects  should  not  be  of  impor¬ 
tance  for  the  degree  of  prediction  accuracy  being  considered.  Finally,  since  the  same 
prediction  program  and  theory  are  used  to  predict  all  nine  satellite  inclinations,  an  error 
in  theory  or  program  should  be  evident  in  all  nine  comparison  plots.  Thus,  it  appears 
that  an  anomaly  exists  in  the  orbital  computations  made  from  the  observations  for  these 
three  satellites  during  1966  and  early  1967. 


PHYSICAL  INTERPRETATION  OF  RESULTS 

The  various  perturbations  in  inclination  indicated  by  the  nine  comparison  plots  can 
be  easily  associated  with  a  particular  change  in  one  or  more  orbital  elements  of  the  sun 
or  moon.  For  example,  an  oscillatory  change  in  inclination  having  a  period  of  about  6 
months  is  evident  in  each  plot.  This  is  related  to  the  apparent  motion  of  the  sun  about 
the  earth  in  1  year  (the  orbital  element  changing  is  the  mean  anomaly  of  the  sun).  It 
should  be  noted  that  the  orbital  elements  of  the  sun  or  moon  appear  as  multiple-frequency 
terms  in  the  inclination  prediction  equation,  hence,  the  oscillation  frequencies  evident  in 
the  graphs  will  be  integral  multiples  of  the  frequencies  of  particular  sun,  moon,  or  satel¬ 
lite  elements. 
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Fig.  6  -  Theoretical  and  observed  inclination  for  Satellite  671 
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Fig.  7  -  Theoretical  and  observed  Inclination  for  Satellite  704 
•  Data  o  Theoretical 


Fig.  8  -  Theoretical  and  observed  inclination  for  Satellite  1514 
•  Data  o  Theoretical 
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Fig.  10  -  Theoretical  and  observed  inclination  for  Satellite  509.  •  Data  o  Theoretical 
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Fig.  11  -  Theoretical  and  observed  inclination  for  Satellite  801.  •  Data  o  Theoretical 
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Fig.  12  -  Theoretical  and  observed  inclination  for  Satellite  1314,  •  Data  o  Theoretical 
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Fig.  13  -  Theoretical  and  observed  inclination  for  Satellite  1420.  *  Data  o  Theoretical 


A  daily  plot  of  predicted  inclination  for  Satellite  902  is  shown  in  Fig.  14.  An  oscilla¬ 
tion  having  a  14 -day  period  is  evident  in  the  graph.  This  is  associated  with  the  motion 
(mean  anomaly  is  the  element  changing)  of  the  moon  around  the  earth  in  28  days.  This 
oscillation  accounts  for  the  apparent  noiselike,  short-period  variations  in  inclination  on 
the  3-year  plots.  Figures  15.  16,  and  17  show  long-term  predictions  in  Inclination  for 
Satellite  902.  Two  major  oscillatory  components  having  periods  of  56  years  and  26  years 
are  evident  in  Fig.  17.  The  first  component  is  due  to  the  rotation  of  the  moon's  orbital 
plane  in  18.6  years,  which  is  the  change  in  right  ascension.  The  second  is  related  to  the 
combined  motions  of  lunar  perigee  and  the  satellite  right  ascension. 


CONCLUSIONS 

The  theoretical  expression  used  in  this  report  for  predicting  the  inclination  of  a 
satellite  orbit  appears  to  yield  fairly  accurate  results  over  a  long  period  of  time.  The 
predicted  data  are  believed  to  be  accurate  to  +0.002 c  for  these  particular  satellites.  To 
obtain  predictions  accurate  to  four  decimal  places  in  inclination,  it  is  necessary  to  increase 
the  complexity  of  the  predicting  program.  Higher  order  terms  that  were  neglected  in  the 
development  of  the  basic  predicting  equation  (Eq.  '1))  must  now  be  included.  For  those 
satellites  having  a  high  surface-area-to-mass  ratio,  it  may  be  necessary  to  correct  for 
solar  radiation  pressure.  The  effects  of  solar  tide  also  become  important  for  this  degree 
of  accuracy  (see  Ref.  7). 

Inclination  is  one  of  the  easiest  satellite  parameters  to  predict.  Prediction  of  right 
ascension  and  argument  of  perigee  are  two  of  the  most  difficult.  In  addition  to  the  effects 
of  the  sun  and  moon  on  these  parameters  as  with  inclination,  the  effects  of  the  earth’s 
gravitational  irregularities  complicate  the  situation.  Though,  as  the  satellite  altitude 
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becomes  high,  these  two  variables  become  easier  to  predict.  Long-term  prediction 
accuracy  for  these  two  parameters  will  increase  as  greater  knowledge  of  the  earth’s 
gravitational  structure  is  obtained. 

SYMBOLS 

a  Semimajor  axis  of  a  satellite  orbit 

ad  Semimajor  axis  of  lunar  or  solar  orbit 

b  si  or  (a  -  nd) 

d  Time  in  ephemeris  days 

e  Eccentricity  of  a  satellite  orbit 

ed  Eccentricity  of  lunar  or  solar  orbit 

Eccentricity  of  the  lunar  orbit 
es  Eccentricity  of  the  sun 

G  Universal  gravitational  constant 

ls  Obliquity  of  the  ecliptic 

i  Inclination  of  a  satellite  orbit 

id  Inclination  of  the  lunar  or  solar  orbit 

Jn  nth  zonal  harmonic  of  the  earth 

Pd 

L  3 

d 

L  Mean  longitude  of  the  sun 

(L  -  I’)  Mean  anomaly  of  the  sun 
u  Mean  anomaly  of  a  satellite 

Mean  anomaly  of  the  sun  or  moon 
NL  Nutation  in  longitude 

n  Mean  motion  of  a  satellite 

p  Total  precession  since  a  given  epoch 

T  Time  in  Julian  centuries  (1  Julian  century  = 

36,525  ephemeris  days) 

a  Inclination  of  the  lunar  orbit  to  the  ecliptic 

Mean  longitude  of  perigee  of  the  sun 


r 
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I'*  Mean  longitude  of  lunar  perigee  measured  along  the 
ecliptic  and  then  along  the  lunar  orbit 
Ms  Nutation  in  obliquity 
<  Mean  obliquity  (inclination)  of  the  ecliptic 
nd  Product  of  G  and  the  mass  of  the  sun  or  moon 
u  Argument  of  perigee  of  a  satellite  orbit 
<od  Argument  of  perigee  of  lunar  or  solar  orbit 
<o  Rate  of  change  of  argument  of  perigee 
fl  Right  ascension  of  the  ascending  node  of  a  satellite  orbit 
nrf  Right  ascension  of  the  sun  or  moon 
n  Rate  of  change  of  right  ascension 
€  Mean  longitude  of  the  moon  measured  along  the 
ecliptic  and  then  along  the  lunar  orbit 
Q  Right  ascension  of  the  lunar  orbit  measured 
along  the  ecliptic 
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Appendix 

FORTRAN  PROGRAM  FOR  PREDICTION  CALCULATIONS 

FTN5.4A  02/24/69 


PR0GRAM  INCkATE 

DIMENSION  DID1  U0U0),D1  DTK  10Qfl),DlDT2(  1000)  »  SAT  1(1000) 

CiDfcLTA(lOOO) 

PR  =  3, 141^926/130, 

T  *  24106,5/36525.  -  1./36525, 

0  *  24106 ,  5  -  1. 

VS1  *  S1N(5.145*DR) 

VCi  =  C0S(5,l45*DR) 

DNT  =  0. 

PK  *  5<l ,  25/36011 ,  *0R 

KM  s  0 

KG  :  II 
DELE  =  0. 

ME  *  0 
VA  =  0 
L  s  li 
IK  c  II 

SEMI  =  1,2027 
W  s  242 , 541.DR 
8  *  21b, 570. DR 
P  =  ,00279 

A  I  =  yt) ,  20<}i>*UI< 

AN  *  12.yl4.360, 

1)0  1»  KO  =  1,10,1 
00  II’-  K  s  l,lfc26-i,l 
T  =  T  ♦  1,/io 525, 

H£  !  )t  *  1 
I)  s  B  ♦  1. 

A  IDs  (23, 452294 -U,  013012s.  T-0, 000150.(  7.. 2  > -0, 01)0003.  (  T»*3  )  )•£ 
ED=((,.Olh7slO4-O,OUOO4l0O.T-O,OOOOOU126.(T*.if ) ) 

WD*(2tl. 22 Ob 3+0. 0000470664»D+0,0004S3»{T**2).U, 000003.(7. .3)) 
WQ  *  Kl)  ♦  ON T 

AMD* (35b. 47563*0 .98560026/ *D-0 ,000150* ( T»*2 ) -0 . UO0003* ( T* *3  ) ) 
RDU  *  1 , 496fc Oil 
AK  *  .9/142 

SLUNG  *  <27 V.bVobb  ♦  0 ,98s6473354*D  ♦  0 . 000303»T»*2) *DR 
SL0NG  =  SLUNG  ♦  DNT 
L  *  l  ♦  1 
«W  *  1*  *  1 

6S  *  ti 

VW  *  0, 

1F<365*KG>40,41,40 

41  KM  *  0 

40  KM  *  KM  +  1 

KG  *  Kl>  ♦  1 

P  s  PK.KM/305, 

FA  =  0, 

f  =  1. 

BP  =  b 

IF(H-3, 14159)43,43,42 

42  HP  s  2, .3, 141s926-b 
FA  *  1. 

43  CONTINUE 

CTH  =  C0S(Alo)*w0S<Al  )*SI.\I(AID1*SIN(A1  >»C0b(UP> 

STH  s  SIN(AC«aSUTM>  ) 

SIR  *  SIM(AJ  )*SINHP)/SIH 
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R  »  AS  IN (SI R ) 

R  *  R  -  FA*2,*P 

EPN  *  AC0S(C0S<AlD*DELe)*CTH-F*SIN<AlD*DELfc>,STH*C0S<»*P)) 

OELTA(L)  =  ( tPN" A  I )  /DR 
AID  a  All)  ♦  DELS 
VS2  s  S I N  ( A 1 1) ) 

VC2  *  CDS (AID) 

2  C0NTJNUE 
HMD  a  aID/2, 

51  =  S1N(HAI0)**4 

52  s  CDS (HA  1 0 ) «*4 
AC  =  7*b'U/2 . 

Al  a  2,»(B  ♦  HD  ♦  AMD) 

A2  *  2,*<B*wD*aMU/2, ) 

43  a  2*b  ♦  2*taD  *3»AMU 
HA  I  a  Al/2, 

A£4  =  2.«d  ♦  2,*WD  *  4,*AMP 

0Rl  *  <l,-5,/2.*fcD**2)*SI\J<Al>-ED/2,*SI'I(A2)*AC*SINU3)*l7,/2.*ED 
C  **2*SlNCAfc4) 

A4  =  2*(B-riD-AMU) 

AS  «  2* ( B-WU-AMD/2 ) 

A6  =  2»b-2«wD-3* AMD 

AE6  a  2.*b  -  2,*RD  -  4 , • AMD 

BRH  a  (l.-5./2.*EU**2)«SlN(A4)-En/2,»SI'I(A5)*AC*SlN(A6>tl7,/2.*ED 
C  **2*S I N ( At6 ) 

P1=SIM  A  I )  * ( Sl«bRi  ♦  S2*BR2  *,5*  ( SlN ( A  I D) **2 ) *SlN< 2 « *B ) » ( 1 . *3 , * 

C  ED»C6S(AMD  )  ♦  3,/2,*fcD**2  ♦  9  , /?.,  *ED**2*C0S  (  2,  *AMD ) ) ) 

P2  a  0, 

T1  a  <l.*o./2,*t**2-2.*E*C0S<AMD)-.!>*E**2*C0S<2.*AMD> )*P1*5,/4.*E 
C  «*2*P2 

U  *  T1  ♦  C9S(Al )»SIN(AID)«( < S !N(HA l D ) ** 2) *S IN( B*2, *WD*2, • AMD) 

C  -(C0S(HA1D)**2)*SIN(8-2,*WD»2,*AMD)  ♦  C0S< AlU)*i>IN(B) ) 

TF  a  (3,»AK/(4,*AN*(l.»t»*2)**i5))*Tl 
G0  TG  (,4,5),VW  ♦  1 
4  C0NTINUE 
VW  a  VW  ♦  1 
DIDTI(L)  a  TF 
MODIFY  M00N  PARAMETERS 
ED  a  0 ,0^4900409 

WD  a  {334.32955o*O,lU4O4O0O3*n-O,OlO35>5*T**2-U,UOQO12*T**3)*DR 
WD  a  wl)  ♦  DNT 

ZL  a  (270.434164*13,l/63965268*D-0,001133*T**2*0f 00Q0019*T**3)*DR 
ZL  a  ZL  +  DNT 

RAS*< 259, 183275-0, 0529539222*0*0. 092078* T* *2* 0,0V0Q02*T* *3 )*DR 
RAS  a  RAS  ♦  DNT 

DNTs(- (17.234*0. 017*T)*SIN(RAS)*0.209*S I N(2*RAS)'1,272*SIN(2*SL0NG 
C  )  -U.204*SIM2*2l>>/36QO,*DR 

DELE  *(<9,210*0, 00  09  •T)*C0S(RAS)-(',!J9O*C3S(2«KAS)*,551*C0S(2* 

C  SL0NG)*O,O88*C0S(2*ZL> )/3600,.dr 

SRAS  =  SIN(RAS) 

CRAS  a  C0S (RAS ) 

AID  a  AC0S ( VC2* VC1  -  VS2*VSl*CRAS> 

SIAJD  a  SIN(AIO) 

RTAS  a  ASIN(VS1*SRAS/SIA10) 

8  *  H  -  HTAS 
GAM  a  WO 
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WD  =  ASH(VS2«SRAS/SIA!I)>  ♦  GAM  -  RAS 
AHU  =  ZL  -  GAH 
AK  *  2.1093 
RDU  =  3.844fc'J0« 

GS  rs  2 

5  C0MTSNUt 

D I D  T  2  <"  U )  =  TF 

D!DT(u  -  Diautu  +  uiuT2a> 

C  M0DIFT  SAT  PARAMETERS 

HOOT  s  ( (6,263l8»0,00l623b*<2,-5.»SlR{A!)**2/2. > J / ( SEM 1**2* 

C  )»AN/>i60, 

OD«r  =  {6.2<j3lri.3,*,OOO:,4115.C0S(Al  )/«S£KI**2*U»-F**2)*.2)  ) 
C  *AN/«J&0, 

BOOT  =  -9D0T 
A I  s  A I  +  01UT(L>*DR 

satkl)  =  ai/,m< 
w  a  *  ♦  Wily  r 
H  =  PS  ♦  BOoT 
609  C0MT INUfc 
100  C0NT  I  (xllt 

IF ( 1000-L ) 8, d , 10 

9  C0>MII,\Ufc 

PRINT  6, (SAT1 (M) ,m=1,10U0#7) 

PRINT  6, (DfcLTA(M) .Msl, 1000,7) 

6  F0Rf1Ar<lX,lO<rx,Fll,6>> 

L  =  0 

10  C0NIUUB 
ifNU 
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This  report  describes  an  investigation  of  the  ability  to  predict  long-term  orbital 
elements  for  medium-height  (satellite  orbiting  the  earth  at  altitudes  above  the 
regions  of  appreciable  drag),  earth-orbiting  satellites.  The  particular  types  of 
satellites  studied  are  those  having  nearly  circular,  polar  orbits.  The  primary  ele¬ 
ment  of  interest  is  orbital  inclination. 

A  theoretical  equation  for  predicting  the  change  in  inclination  caused  by  the 
gravitational  attraction  of  the  sun  and  moon  is  given.  Expressions  for  the  orbital 
elements  of  the  sun  and  moon  required  in  this  inclination- predicting  equation  are 
developed.  A  FORTRAN  program  listing  and  computation  procedure  for  determining 
the  inclination  of  a  satellite  at  any  time  after  a  chosen  epoch  is  in  the  appendix. 

Nine  satellites  were  selected  for  study.  For  these  satellites,  comparison  plots 
of  observed  and  predicted  inclinations  over  a  3-year  time  interval  are  shown.  The 
agreement  between  predicted  and  observed  inclinations  is  very  good  for  most  of  the 
satellites.  Theoretical  plots  of  the  inclination  of  Satellite  902  over  time  intervals  of 
u/>  to  300  years  are  included. 
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